ptm = proc.time()Frequentemente utilizamos serviços de classificação, seja com o intuito de recomendar/sugerir algo a um usuário final, como Netflix, Mercado Livre, OLX dentre outros, ou seja para determinar em testes de diagnóstico de doenças se um paciente deve ser classificado como doente, não doente ou suspeito. Além de outros exemplos. O conjunto de dados explanados a seguir foram extraídos da plataforma Kaggle e são referentes ao catálogo de filmes da plataforma IMDB, o intuito aqui é criar um sistema de classificação de estrelas do catálogo de avaliação de filmes do IMDB.
Esse projeto irá explorar o conjunto de dados a fim de gerar insights, e aplicar o algoritmo de Random Forest - RF (baseado na clusterização de Breiman) para classificação de filmes nas seguintes faixas/categorias de estrelas (Rating):
Para este, se faz necessário um largo conjunto de dados que, por sua vez, deve ser bem estruturado e consistente. Portanto, é fundamental manipular, limpar e remover informações faltantes dos dados para a implementação do modelo de Machine Learning do Random Forest que será treinado com orientação a esses dados.
Além disso, esse trabalho irá levantar quais os fatores mais importantes para que um filme tenha uma alta classificação (categoria de faixas de estrelas), via modelagem de Random Forest, de tal forma que essas mesmas variáveis/fatores destacados como importantes no algoritmo do RF possam ser aproveitados posteriormente em análises futuras de classificação. Resultados e script foram gerados utilizando recursos e linguagens em \(R\), \(\LaTeX\) e \(Markdown\).
Pacotes e Leitura
Alguns pacotes R utilizados nessa rotina
library(tidyverse) # Manipulação de banco de dados e análise exploratória
library(highcharter) # visualização Gráfica
library(DT) # Construção de Tabelas
library(knitr) # Opções Rmarkdown (inclusão de tabelas, imagens, etc.)
library(kableExtra) # Construção de Tabelas
library(rpart.plot) # Recursive Partitioning and Regression Trees
library(randomForest) # Modelo Random Forest
library(corrplot) # Matriz de correlação - Visualização gráfica
library(wordcloud2) # Numvem de palavras (Word Cloud)
library(caret) # Matriz de confusão (Confusion Matrix)
library(rpart)
#devtools::install_github("ropenscilabs/icon")
library(icon)RECOMEND.METADATA = readxl::read_xlsx(path = ".../DB/IMDB.xlsx", sheet = 1)
db = RECOMEND.METADATA+ + =
O conjunto de dados utilizado foi o de título IMDB 5000 extraído da plataforma kaggle. As informações contidas no banco foram catalogadas de filmes publicados ao longo de 100 anos em 66 países (entre 1916 e 2016) da plataforma IMDB - Internet Movie DataBases, o Arquivo original contém 5044 filmes (observações) e 28 variáveis descritas a seguir.
Os dados em questão são públicos e disponíveis para download clicando AQUI.
## Warning in instance$preRenderHook(instance): It seems your data is too
## big for client-side DataTables. You may consider server-side processing:
## https://rstudio.github.io/DT/server.html
## O banco de dados possui 5043 observações e 28 variáveis
Excluindo algumas variáveis que não serão utilizadas.
O motivo de remover as variáveis em (1) é devido à tentativa de reduzir o número de missing na base de dados. no passo 3.1 de Tratamento de NA’s.
#db=RECOMMEND.METADATA
db = db %>% select(- director_facebook_likes, -actor_3_facebook_likes, -actor_2_name, -actor_1_facebook_likes, -actor_1_name, -actor_3_name,-facenumber_in_poster, -actor_2_facebook_likes, -color, -director_name, -language, -country, -content_rating, -movie_imdb_link)## O banco de dados agora possui 14 variáveis. Ou seja, foram removidas 14 variáveis no passo anterior.
Estrutura dos dados
glimpse(db) #str(db)## Observations: 5,043
## Variables: 14
## $ num_critic_for_reviews <int> 723, 302, 602, 813, NA, 462, 392, 32...
## $ duration <int> 178, 169, 148, 164, NA, 132, 156, 10...
## $ gross <int> 760505847, 309404152, 200074175, 448...
## $ genres <fct> Action|Adventure|Fantasy|Sci-Fi, Act...
## $ movie_title <fct> Avatar , Pirates of the Caribbean: A...
## $ num_voted_users <int> 886204, 471220, 275868, 1144337, 8, ...
## $ cast_total_facebook_likes <int> 4834, 48350, 11700, 106759, 143, 187...
## $ plot_keywords <fct> avatar|future|marine|native|parapleg...
## $ num_user_for_reviews <int> 3054, 1238, 994, 2701, NA, 738, 1902...
## $ budget <dbl> 237000000, 300000000, 245000000, 250...
## $ title_year <int> 2009, 2007, 2015, 2012, NA, 2012, 20...
## $ imdb_score <fct> 7.9, 7.1, 6.8, 8.5, 7.1, 6.6, 6.2, 7...
## $ aspect_ratio <fct> 1.78, 2.35, 2.35, 2.35, , 2.35, 2.35...
## $ movie_facebook_likes <int> 33000, 0, 85000, 164000, 0, 24000, 0...
Transformando as variáveis imdb_score e aspect_ratio em numéricas
db$imdb_score = as.numeric(as.character(db$imdb_score))
db$aspect_ratio = as.numeric(as.character(db$aspect_ratio))hc1 = hchart(db$imdb_score, color = "#e8bb0b", name = "imdb_score") %>%
hc_title(text = "Histograma dos votos na plataforma IMDB") %>%
hc_exporting(enabled = TRUE, filename = "Fig1-Pimenta"); hc1hc2 = hchart(db$num_voted_users, color = "#786eea", name = "imdb_num_votos") %>%
hc_title(text = "Histograma dos número de votos na plataforma IMDB") %>%
hc_exporting(enabled = TRUE, filename = "Fig2-Pimenta"); hc2Criando uma função para cálculo da moda
Moda <- function(x) {
ux <- unique(x)
ux[which.max(tabulate(match(x, ux)))]
}## Média, moda e mediana da variável imdb_score são, respectivamente, 6.44 6.7 6.6
Extraindo valores da variável gênero e transformando em dummies
gg <- as.character(db$genres); gg <- gsub("-", "_", as.character(gg))
#t <- unlist(strsplit(gg[1],split = "\\|"))
tem1 <- data.frame()
for(i in 1:length(gg)){
tem <- tem1
t <- unlist(strsplit(gg[i],split = "\\|"))
temp <- data.frame(t)
tem1 <- rbind(tem,temp)
}
Gen <- unique(tem1); Gen <- gsub("-", "_", as.character(Gen$t))
cat("Existem", length(Gen), "valores de gêneros únicos de filme no banco de dados.")## Existem 26 valores de gêneros únicos de filme no banco de dados.
Genname <- Gen
fe <- matrix(data = 0, nrow = length(gg), ncol = length(Genname))
fe <- data.frame(fe); colnames(fe) <- Genname
i=j=0
for(i in 1:length(gg)){
for(j in 1:length(Genname)){
g <- grepl(Genname[j], gg[i])
if(g == TRUE){
fe[i, j] <- 1
}
}
}
NumGen = as_tibble(rbind(apply(fe,2,sum)))
NumGen = gather(NumGen, key = "variables", value = "num_gender")
NumGen = NumGen[order(NumGen$num_gender, decreasing = TRUE), ]
hc3 <- highchart() %>%
hc_add_series(data = NumGen$num_gender,
type = "bar",
name = "# de filmes",
showInLegend = FALSE,
tooltip = list(valueDecimals = 0, valuePrefix = "", valueSuffix = ""), color="blue") %>%
hc_yAxis(title = list(text = "Quantitativo de filmes"),
allowDecimals = TRUE, max = (max(NumGen$num_gender)+103),
labels = list(format = "{value}")) %>%
hc_xAxis(title = list(text = "Gênero de filme"),
categories = NumGen$variables,
tickmarkPlacement = "on",
opposite = FALSE) %>%
hc_title(text = "Quantitativo de filmes por gênero",
style = list(fontWeight = "bold")) %>%
hc_subtitle(text = paste("")) %>%
hc_tooltip(valueDecimals = 2,
pointFormat = "{point.y} filmes")%>%
#pointFormat = "Variável: {point.x} <br> Missing: {point.y}")
hc_credits(enabled = TRUE,
text = "Fonte: IMDB/KAGGLE. Elaboração: Ewerson Pimenta.",
style = list(fontSize = "10px")) %>%
hc_exporting(enabled = TRUE, filename = "F3-filmes-genero-Pimenta")
#hc <- hc %>%
# hc_add_theme(hc_theme_darkunica())
hc3; remove(gg,t,tem, temp, tem1,Gen, g, i, j)Algumas variáveis de gêneros possuem pouquíssimos filmes serão removidas do banco. São elas Film_Noir, Short, News, Reality_TV e Game_Show.
fe$Film_Noir <- NULL
fe$Short <- NULL
fe$News <- NULL
fe$Reality_TV <- NULL
fe$Game_Show <- NULLUnificando a base de dados com as novas variáveis de gêneros únicos descobertos.
db1 = as_tibble(cbind(as.data.frame(db),fe))
cat("Foram inseridas ", dim(fe)[2],"novas variáveis provenientes dos gêneros únicos descobertos no passo anterior.")## Foram inseridas 21 novas variáveis provenientes dos gêneros únicos descobertos no passo anterior.
Semelhante à variável de gênero, foi feito o split por palavra chave da variável keyword, esse, por sua vez, continha 6780 palavras únicas. E então, geramos uma nuvem de palavras chave, com o auxílio do pacote wordcloud2 com o seguinte comando, após obter o data frame do nome das palavras chaves e frequência.
library(wordcloud2)
wordcloud2(NumKeyWord)Excluindo filmes lançados antes de 1980
hc00 = hchart(db$title_year, color = "#53a074", name = "imdb1_score") %>%
hc_title(text = "Histograma do ano de publicação dos filmes") %>%
hc_exporting(enabled = TRUE, filename = "Fig1-Pimenta"); hc00Percebemos pelo gráfico anterior que existem poucos filmes publicados antes de 1980 e estes podem não ser representativos. Decidimos, então, por trabalhar apenas com os filmes que foram publicados a aprtir de 1980
db <- db[db$title_year >= 1980,]
hc01 = hchart(db$title_year, color = "#79d8a2", name = "imdb1_score") %>%
hc_title(text = "Histograma do ano de publicação dos filmes") %>%
hc_exporting(enabled = TRUE, filename = "Fig1-Pimenta"); hc01Porcentagem de NA por variável
db_miss <- db %>% summarise_all(funs(sum(is.na(.))/n()))## Warning: package 'bindrcpp' was built under R version 3.5.1
db_miss <- gather(db_miss, key = "variables", value = "percent_missing")
db_miss$percent_missing = 100*db_miss$percent_missing
db_miss = db_miss[order(db_miss$percent_missing, decreasing = TRUE), ]
#db_miss
hc4 <- highchart() %>%
hc_add_series(data = db_miss$percent_missing,
type = "bar",
name = "Porcentagem de missing",
showInLegend = FALSE,
tooltip = list(valueDecimals = 2, valuePrefix = "", valueSuffix = " %")) %>%
hc_yAxis(title = list(text = "Porcentagem de missing"),
allowDecimals = TRUE, max = 20,
labels = list(format = "{value}%")) %>%
hc_xAxis(categories = db_miss$variables,
tickmarkPlacement = "on",
opposite = FALSE) %>%
hc_title(text = "Porcentagem de missinig por variável",
style = list(fontWeight = "bold")) %>%
hc_subtitle(text = paste("")) %>%
hc_tooltip(valueDecimals = 2,
pointFormat = "Missing: {point.y}")%>%
#pointFormat = "Variável: {point.x} <br> Missing: {point.y}")
hc_credits(enabled = TRUE,
text = "Fonte: IMDB/KAGGLE. Elaboração: Ewerson Pimenta.",
style = list(fontSize = "10px")) %>%
hc_exporting(enabled = TRUE, filename = "Fig0-Pimenta")
#hc <- hc %>%
# hc_add_theme(hc_theme_darkunica())
hc4Removendo todas as observações que contêm NA.
db1 <- db1 %>% drop_na(gross, budget, aspect_ratio
,title_year
,num_critic_for_reviews, num_user_for_reviews)
#glimpse(db1)
db1 <- db1 %>% na.omit()
db1 <- db1 %>% drop_na()Removendo dados duplicados
cat("Existem", sum(duplicated(db1)), "filmes duplicados na base de dados.")## Existem 33 filmes duplicados na base de dados.
db1 <- db1[!duplicated(db1), ]## O novo banco de dados, sem observações faltantes, possui 3782 observações. Ou seja, no processo de remoção de valores faltantes foram perdidas 1000 observações. Finalmente, removemos todas as observações duplicadas e faltantes do banco.
Porcentagem de ZEROS por variável
zeros <- (colSums(db1==0)/nrow(db1)*100); var <- names(zeros)
db_zero <- data.frame(var,zeros); rownames(db_zero) <- NULL
db_zero <- db_zero[order(db_zero$zeros, decreasing = TRUE), ]
hc4_1 <- highchart() %>%
hc_add_series(data = db_zero$zeros,
type = "bar",
name = "Porcentagem de zeros",
showInLegend = FALSE,
tooltip = list(valueDecimals = 2, valuePrefix = "", valueSuffix = " %"), color="pink") %>%
hc_yAxis(title = list(text = "Porcentagem de zero"),
allowDecimals = TRUE, max = 100,
labels = list(format = "{value}%")) %>%
hc_xAxis(categories = db_zero$var,
tickmarkPlacement = "on",
opposite = FALSE) %>%
hc_title(text = "Porcentagem de zeros por variável",
style = list(fontWeight = "bold")) %>%
hc_subtitle(text = paste("")) %>%
hc_tooltip(valueDecimals = 2,
pointFormat = "Zeros: {point.y}")%>%
#pointFormat = "Variável: {point.x} <br> Missing: {point.y}")
hc_credits(enabled = TRUE,
text = "Fonte: IMDB/KAGGLE. Elaboração: Ewerson Pimenta.",
style = list(fontSize = "10px")) %>%
hc_exporting(enabled = TRUE, filename = "Fig00-Pimenta")
#hc <- hc %>%
# hc_add_theme(hc_theme_darkunica())
hc4_1Aqui podemos perceber que os zeros estão concentrados nas variáveis de gênero que são binárias e não faz sentido analizar os zeros nessas variáveis.
paste0("Entretanto, pode-se perceber uma proporção de zeros na variável movie_facebook_likes igual a ",round(db_zero[22,2],2), "% e, também na variável cast_total_facebook_likes de ",round(db_zero[23,2],2),"%. Decidimos por remover a variável movie_facebook_likes.")## [1] "Entretanto, pode-se perceber uma proporção de zeros na variável movie_facebook_likes igual a 46.17% e, também na variável cast_total_facebook_likes de 0.13%. Decidimos por remover a variável movie_facebook_likes."
db1$movie_facebook_likes <- NULLUm passo importante é penalizar a variável de escore imdb_score pelo número de votos recebidos. Ver mais no estimador de Shrinkage
\(WR = \frac{v}{v+m} \times R + \frac{m}{v+m} \times C\)
Onde,
R = as.numeric(db1$imdb_score)
v = as.numeric(db1$num_voted_users)
m = summary(db1$num_voted_users)[2] # 1st Qu.
C = mean(db1$imdb_score)
db1$WR = (v/(v+m))*R + (m/(v+m))*CMedidas pontuais e de dispersão de imdb_score e WR
summary(db1$imdb_score)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 1.600 5.900 6.600 6.468 7.200 9.300
summary(db1$WR)## Min. 1st Qu. Median Mean 3rd Qu. Max.
## 2.578 6.166 6.498 6.573 6.983 9.269
cat("Além disso, conseguimos reduzir o desvio padrão que orignalmente era de ",round(sd(db1$imdb_score),2),",para a variável imdb_score, e agora passou a ser ",round(sd(db1$WR),2), "para a variável WR.")## Além disso, conseguimos reduzir o desvio padrão que orignalmente era de 1.05 ,para a variável imdb_score, e agora passou a ser 0.72 para a variável WR.
Selecionando uma amostra aleatória (a.a.) de tamanho \(n=800\) e representando em um gráfico de dispersão com valores reais vs ajustados.
set.seed(123654)
amostra0 = sample(x = 1:dim(db1)[1], size = 800, replace = FALSE)
dbX = db1[amostra0,] %>%
select(movie_title,num_voted_users,imdb_score, WR)
dss <- map(c("cross"), function(s){
x <- as.numeric(dbX$imdb_score)
y <- as.numeric(dbX$WR)
list(name = s,
data = list_parse(data_frame(x, y)),
marker = list(symbol = s, enabled = TRUE), lineColor = "#56667a")
})
#dss[[1]]$data[amostra1]
hc5 = highchart() %>%
hc_chart(type = "scatter", color = "#56667a") %>%
hc_title(text = "Score IMDB vs WR (calculado pelo estimador de Shrinkage)") %>%
hc_subtitle(text = "800 filmes selecionados via amostra aleatória simples") %>%
hc_xAxis(title = list(text = "x: imdb_score"),
allowDecimals = TRUE, labels = list(format = "{value}★")) %>%
hc_yAxis(title = list(text = "y: WR (calibrado)"),
allowDecimals = TRUE, labels = list(format = "{value}★")) %>%
hc_exporting(enabled = TRUE, filename = "F3-Pimenta") %>%
hc_add_series_list(dss); hc5; remove(dbX, dss)Análise do cálculo de IMDB
Ou seja, para os filmes catalogados com m muito inferior a 18mil o novo escore calibrado teve maior diferentça que aqueles superior a 18 mil. Além disso, quando scores são maiores que 6,5 o WR tende a cair, caso contrário o valor pode aumentar. Quanto maior o número de votos recebidos, menor a diferença do valor de IMDB_score e WR.
set.seed(123654)
amostra2 = sample(x = 1:dim(db1)[1], size = 35, replace = FALSE)
db1[amostra2,] %>%
select(movie_title,num_voted_users,imdb_score, WR, budget) %>%
datatable(rownames = amostra2,options = list(searchin = TRUE, scrollX = TRUE, pageLength = 5))Finalmente, segue uma a.a.da base de dados após limpeza e análise exploratória.
set.seed(123851)
amostra3 = sample(x = 1:dim(db1)[1], size = 35, replace = FALSE)
db1[amostra3,] %>%
datatable(rownames = amostra3,options = list(searchin = TRUE, scrollX = TRUE, pageLength = 5))Para os dados quantitativos a seguinte Matriz de correlação. Foi gerada utilizando o pacote corrplot
#install.packages("corrplot")
library(corrplot)
res1 <- cor.mtest(db_cor, conf.level = .95)
res2 <- cor.mtest(db_cor, conf.level = .99)
corrplot::corrplot(cor(db_cor), method = "color", col = col(200),
type = "upper", order = "hclust", number.cex = .7,
addCoef.col = "black", # Add coefficient of correlation
tl.col = "black", tl.srt = 90, # Text label color and rotation
# Combine with significance
p.mat = res1$p, sig.level = 0.01, insig = "blank",
# hide correlation coefficient on the principal diagonal
diag = FALSE)Descrição 1. Random Forest foi desenvolvido para agregar árvores de decisão (modelo de classificação);
2. Pode ser usado para modelo de classificação (p/ var. resposta categórica) ou regressão (no caso de haver variável resposta contínua);
3. Evita overfitting;
4. Permite trabalhar com um largo número de características de um conjunto de dados;
5. Auxilia na seleção de variáveis baseada em um algoritmo que calcula a importância por variável (assim, tendo conhecimento de quais variáveis são mais importantes, podemos usar essa informação para outros modelos de classificação);
6. User-friendly: apenas 2 parâmetros livres:
Passo-a-Passo
É realizado em 3 passos:
Exemplo
Transformação da variável contínua WRem categórica WR_Grp
Transformando a variável WR em fator, com as seguintes categorias:
movie = db1
Grp <- function(tn){
tn = abs(tn)
if (tn >= 0 & tn < 4){
return('[0, 4)')
}else if(tn >= 4 & tn < 6){
return('[4, 6)')
}else if(tn >= 6 & tn < 8){
return('[6, 8)')
}else if (tn >= 8){
return('de 8 a 10')
}
}
# apply the Group function to the WR column
movie$WR_Grp <- sapply(movie$WR,Grp)
# set as factor the new column
movie$WR_Grp <- as.factor(movie$WR_Grp)
#View(head(movie))
table(movie$WR_Grp)##
## [0, 4) [4, 6) [6, 8) de 8 a 10
## 13 661 2974 134
# apply the Group function to the WR column
imdb_score_Grp <- sapply(movie$imdb_score,Grp)
# set as factor the new column
imdb_score_Grp <- as.factor(imdb_score_Grp)
#View(head(movie))
table(imdb_score_Grp)## imdb_score_Grp
## [0, 4) [4, 6) [6, 8) de 8 a 10
## 84 964 2521 213
Remoção de variáveis Remove as variáveis imdb_score, genres, plot_keywords, movie_imdb_link e WR.
movie$imdb_score <- NULL
movie$genres <- NULL
movie$plot_keywords <- NULL
movie$movie_imdb_link <- NULL
movie$WR <- NULLglimpse(movie)## Observations: 3,782
## Variables: 32
## $ num_critic_for_reviews <int> 723, 302, 602, 813, 462, 392, 324, 6...
## $ duration <int> 178, 169, 148, 164, 132, 156, 100, 1...
## $ gross <int> 760505847, 309404152, 200074175, 448...
## $ movie_title <fct> Avatar , Pirates of the Caribbean: A...
## $ num_voted_users <int> 886204, 471220, 275868, 1144337, 212...
## $ cast_total_facebook_likes <int> 4834, 48350, 11700, 106759, 1873, 46...
## $ num_user_for_reviews <int> 3054, 1238, 994, 2701, 738, 1902, 38...
## $ budget <dbl> 237000000, 300000000, 245000000, 250...
## $ title_year <int> 2009, 2007, 2015, 2012, 2012, 2007, ...
## $ aspect_ratio <dbl> 1.78, 2.35, 2.35, 2.35, 2.35, 2.35, ...
## $ Action <dbl> 1, 1, 1, 1, 1, 1, 0, 1, 0, 1, 1, 1, ...
## $ Adventure <dbl> 1, 1, 1, 0, 1, 1, 1, 1, 1, 1, 1, 1, ...
## $ Fantasy <dbl> 1, 1, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, ...
## $ Sci_Fi <dbl> 1, 0, 0, 0, 1, 0, 0, 1, 0, 1, 1, 0, ...
## $ Thriller <dbl> 0, 0, 1, 1, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ Documentary <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ Romance <dbl> 0, 0, 0, 0, 0, 1, 1, 0, 0, 0, 0, 0, ...
## $ Animation <dbl> 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, ...
## $ Comedy <dbl> 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, ...
## $ Family <dbl> 0, 0, 0, 0, 0, 0, 1, 0, 1, 0, 0, 0, ...
## $ Musical <dbl> 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, ...
## $ Mystery <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, ...
## $ Western <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ Drama <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ History <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ Sport <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ Crime <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ Horror <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ War <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ Biography <dbl> 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, 0, ...
## $ Music <dbl> 0, 0, 0, 0, 0, 0, 1, 0, 0, 0, 0, 0, ...
## $ WR_Grp <fct> [6, 8), [6, 8), [6, 8), de 8 a 10, [...
cat("A base de dados que iremos trabalhar no modelo tem", dim(movie)[1] , "observações e ", dim(movie)[2]-1, "variáveis, sem contar a variável que contém os títulos dos filmes.")## A base de dados que iremos trabalhar no modelo tem 3782 observações e 31 variáveis, sem contar a variável que contém os títulos dos filmes.
Partição dos dados
Particionando a base de dados em Treino e Teste, esses dois (Treino e Teste) também terão armazenos os nomes dos filmes selecionados via amostragem probabilística dos dados originais separadamente das bases de Treino e Teste.
Posteriormente, removemos os rótulos dos filmes nas bases Treino e Teste
set.seed(9182345)
ind <- sample(2, nrow(movie), replace = T, prob = c(0.7, 0.3))
train <- movie[ind==1, -4]
test <- movie[ind==2, -4]
trainMovie <- movie[ind==1, 4]
testMovie <- movie[ind==2, 4]Distribuição dos escores nas bases de treino e teste
round(table(train$WR_Grp)/sum(table(train$WR_Grp))*100,2)##
## [0, 4) [4, 6) [6, 8) de 8 a 10
## 0.37 17.74 78.30 3.59
round(table(test$WR_Grp)/sum(table(test$WR_Grp))*100,2)##
## [0, 4) [4, 6) [6, 8) de 8 a 10
## 0.27 16.83 79.46 3.44
Inicialmente utilizaremos o pacote randomForest que implmenta o algoritmo de Random Forest de Breiman (baseado na clusterização de Breiman, originalmente codificada em Fortran) que tem por finalidade classificar e/ou criar regressão. Além disso, pode ser usado em um modelo não supervisionado para avaliar proximidades entre pontos.
Estamos usando, a partir daqui, a base de treino.
#library(randomForest)
#library(rpart)
#library(rpart.plot)
#rf <- randomForest(proximity = T,ntree = 38,do.trace = T,WR~.,data=training)
set.seed(9984512)
rf <- randomForest(WR_Grp~.,data=train)
rf##
## Call:
## randomForest(formula = WR_Grp ~ ., data = train)
## Type of random forest: classification
## Number of trees: 500
## No. of variables tried at each split: 5
##
## OOB estimate of error rate: 15.43%
## Confusion matrix:
## [0, 4) [4, 6) [6, 8) de 8 a 10 class.error
## [0, 4) 0 4 6 0 1.00000000
## [4, 6) 0 184 291 0 0.61263158
## [6, 8) 0 63 2027 6 0.03291985
## de 8 a 10 0 0 43 53 0.44791667
attributes(rf)## $names
## [1] "call" "type" "predicted"
## [4] "err.rate" "confusion" "votes"
## [7] "oob.times" "classes" "importance"
## [10] "importanceSD" "localImportance" "proximity"
## [13] "ntree" "mtry" "forest"
## [16] "y" "test" "inbag"
## [19] "terms"
##
## $class
## [1] "randomForest.formula" "randomForest"
Olhando as 6 primeiras observações real X predito
p1 <- predict(rf,train)
head(p1)## 1 2 3 4 5 6
## [6, 8) [6, 8) [6, 8) de 8 a 10 [6, 8) [6, 8)
## Levels: [0, 4) [4, 6) [6, 8) de 8 a 10
head(train$WR_Grp)## [1] [6, 8) [6, 8) [6, 8) de 8 a 10 [6, 8) [6, 8)
## Levels: [0, 4) [4, 6) [6, 8) de 8 a 10
Matriz de confusão
library(caret)
library(e1071)## Warning: package 'e1071' was built under R version 3.5.1
confusionMatrix(p1, train$WR_Grp)## Confusion Matrix and Statistics
##
## Reference
## Prediction [0, 4) [4, 6) [6, 8) de 8 a 10
## [0, 4) 10 0 0 0
## [4, 6) 0 475 0 0
## [6, 8) 0 0 2096 0
## de 8 a 10 0 0 0 96
##
## Overall Statistics
##
## Accuracy : 1
## 95% CI : (0.9986, 1)
## No Information Rate : 0.783
## P-Value [Acc > NIR] : < 2.2e-16
##
## Kappa : 1
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: [0, 4) Class: [4, 6) Class: [6, 8)
## Sensitivity 1.000000 1.0000 1.000
## Specificity 1.000000 1.0000 1.000
## Pos Pred Value 1.000000 1.0000 1.000
## Neg Pred Value 1.000000 1.0000 1.000
## Prevalence 0.003736 0.1774 0.783
## Detection Rate 0.003736 0.1774 0.783
## Detection Prevalence 0.003736 0.1774 0.783
## Balanced Accuracy 1.000000 1.0000 1.000
## Class: de 8 a 10
## Sensitivity 1.00000
## Specificity 1.00000
## Pos Pred Value 1.00000
## Neg Pred Value 1.00000
## Prevalence 0.03586
## Detection Rate 0.03586
## Detection Prevalence 0.03586
## Balanced Accuracy 1.00000
RF_importance = randomForest::importance(rf)[order(randomForest::importance(rf)[,1], decreasing = TRUE), ]
knitr::kable(RF_importance)| x | |
|---|---|
| num_voted_users | 162.5365055 |
| gross | 88.2600082 |
| num_user_for_reviews | 87.5743810 |
| duration | 87.3916416 |
| budget | 78.6677042 |
| num_critic_for_reviews | 74.3321112 |
| cast_total_facebook_likes | 73.5806642 |
| title_year | 67.2657353 |
| Drama | 34.9327412 |
| Horror | 26.4569613 |
| aspect_ratio | 16.4930280 |
| Action | 14.0246971 |
| Comedy | 12.0400729 |
| Sci_Fi | 11.0373238 |
| Thriller | 10.6637195 |
| Fantasy | 9.7803018 |
| Crime | 9.6097769 |
| Adventure | 9.1179369 |
| Romance | 9.0227241 |
| Animation | 7.6715822 |
| Family | 7.3746822 |
| Mystery | 6.2644388 |
| Music | 4.2544050 |
| Biography | 2.3292811 |
| Sport | 2.2763517 |
| Musical | 1.9174621 |
| War | 1.9038189 |
| Documentary | 1.6206763 |
| Western | 1.5890124 |
| History | 0.8901083 |
| Verificamos que as variávei | Game_Show, Sci_Fi, Reality_TV, News e Film_Noir não foram relevantes para o algoritmo do random forest. |
randomForest::varImpPlot(rf)Taxa de Erro - Random Forest
plot(rf)
legend('topright', colnames(rf$err.rate), col=1:5, fill=1:5)Observamos que, a partir do número de árvores geradas \(ntrees > 60\) o erro OOB (Out of Bag) não pode ser melhorado.
Ajustando e melhorando estimativas
Além disso, iremos alterar alguns parâmetros da função randomForest como o número de ntrees e mtry. Assim, repetimos o algoritmo do Random Forest, ainda usando a base treino.
# Tune mtry
x = as.data.frame(train1[,-31])
y = (as.factor(train1$WR_Grp))
t <- tuneRF(x = x, y = y,
stepFactor = 1.8,
plot = TRUE,
ntreeTry = 60,
trace = TRUE,
improve = 0.05)## mtry = 5 OOB error = 16.1%
## Searching left ...
## mtry = 3 OOB error = 16.96%
## -0.05336427 0.05
## Searching right ...
## mtry = 9 OOB error = 15.84%
## 0.0162413 0.05
Aparentemente, 9 é um bom candidato ao valor de \(m_{try}\)
#set.seed(093180)
set.seed(998451)
rf1 <- randomForest(WR_Grp~.,data=train1,
ntree = 60,
mtry = 9,
importance = TRUE,
proximity = TRUE)
rf1##
## Call:
## randomForest(formula = WR_Grp ~ ., data = train1, ntree = 60, mtry = 9, importance = TRUE, proximity = TRUE)
## Type of random forest: classification
## Number of trees: 60
## No. of variables tried at each split: 9
##
## OOB estimate of error rate: 15.69%
## Confusion matrix:
## [0, 4) [4, 6) [6, 8) de 8 a 10 class.error
## [0, 4) 0 4 6 0 1.00000000
## [4, 6) 0 203 272 0 0.57263158
## [6, 8) 0 92 1998 6 0.04675573
## de 8 a 10 0 0 40 56 0.41666667
RF_importance1 = randomForest::importance(rf1)[order(randomForest::importance(rf1)[,1], decreasing = TRUE), ]
plot(rf1)
legend('topright', colnames(rf$err.rate), col=1:5, fill=1:5)RF1 = data.frame(variables = rownames(RF_importance1), importance = RF_importance1[,6])
RF1 = RF1[order(RF1$importance, decreasing = TRUE),]
rownames(RF1) <- NULL
hc6 <- highchart() %>%
hc_add_series(data = RF1$importance,
type = "bar",
name = "Importância",
showInLegend = FALSE,
tooltip = list(valueDecimals = 2, valuePrefix = "", valueSuffix = "")) %>%
hc_yAxis(title = list(text = "Importância"),
allowDecimals = TRUE, max = 200,
labels = list(format = "{value}")) %>%
hc_xAxis(title = list(text = "Fatores"),
categories = RF1$variables,
tickmarkPlacement = "on",
opposite = FALSE) %>%
hc_title(text = "Importância por fator - Random Forest",
style = list(fontWeight = "bold")) %>%
hc_subtitle(text = paste("")) %>%
hc_tooltip(valueDecimals = 2,
pointFormat = "Importância: {point.y}")%>%
#pointFormat = "Variável: {point.x} <br> Importância: {point.y}")
hc_credits(enabled = TRUE,
text = "Fonte: IMDB/KAGGLE. Elaboração: Ewerson Pimenta.",
style = list(fontSize = "10px")) %>%
hc_exporting(enabled = TRUE, filename = "F4-missing-Pimenta")
#hc <- hc %>%
# hc_add_theme(hc_theme_darkunica())
hc6Removendo variáveis com pouca importância
Decidimos por retirar todos os fatores que retornaram importância abaixo de 10 em MeanDecreaseGini. Portanto, ficaremos apenas com as seguintes variáveis:
RF1 = RF1[1:13,]
hc6_1 <- highchart() %>%
hc_add_series(data = RF1$importance,
type = "bar",
name = "Importância",
showInLegend = FALSE,
tooltip = list(valueDecimals = 2, valuePrefix = "", valueSuffix = "")) %>%
hc_yAxis(title = list(text = "Importância"),
allowDecimals = TRUE, max = 200,
labels = list(format = "{value}")) %>%
hc_xAxis(title = list(text = "Fatores"),
categories = RF1$variables,
tickmarkPlacement = "on",
opposite = FALSE) %>%
hc_title(text = "Importância por fator - Random Forest",
style = list(fontWeight = "bold")) %>%
hc_subtitle(text = paste("")) %>%
hc_tooltip(valueDecimals = 2,
pointFormat = "Importância: {point.y}")%>%
#pointFormat = "Variável: {point.x} <br> Importância: {point.y}")
hc_credits(enabled = TRUE,
text = "Fonte: IMDB/KAGGLE. Elaboração: Ewerson Pimenta.",
style = list(fontSize = "10px")) %>%
hc_exporting(enabled = TRUE, filename = "F6_1-importance-Pimenta")
#hc <- hc %>%
# hc_add_theme(hc_theme_darkunica())
hc6_1Selecionando variáveis importantes
train2 = train %>% select(num_voted_users, gross, duration, num_user_for_reviews, budget, cast_total_facebook_likes, title_year, num_critic_for_reviews, Drama, Horror, aspect_ratio, Action, Comedy, WR_Grp) %>% droplevels()
test2 = test %>% select(num_voted_users, gross, duration, num_user_for_reviews, budget, cast_total_facebook_likes, title_year, num_critic_for_reviews, Drama, Horror, aspect_ratio, Action, Comedy, WR_Grp) %>% droplevels()
movie = rbind(train2, test2)# Tune mtry
x = as.data.frame(train2[,-14])
y = (as.factor(train2$WR_Grp))
t <- tuneRF(x = x, y = y,
stepFactor = .7,
plot = TRUE,
ntreeTry = 60,
trace = TRUE,
improve = 0.05)## mtry = 3 OOB error = 16.55%
## Searching left ...
## mtry = 5 OOB error = 15.95%
## 0.03611738 0.05
## Searching right ...
## mtry = 2 OOB error = 17.03%
## -0.02934537 0.05
Modelo Final Os parâmetros utilizados, finalmente serão \(m_{try} = 5\) e \(n_{tree} = 60\)
set.seed(093180)
rf_final <- randomForest(WR_Grp~.,data = train2,
ntree = 60,
mtry = 5,
importance = TRUE,
proximity = TRUE)
rf_final; #attributes(rf_final)##
## Call:
## randomForest(formula = WR_Grp ~ ., data = train2, ntree = 60, mtry = 5, importance = TRUE, proximity = TRUE)
## Type of random forest: classification
## Number of trees: 60
## No. of variables tried at each split: 5
##
## OOB estimate of error rate: 16.03%
## Confusion matrix:
## [0, 4) [4, 6) [6, 8) de 8 a 10 class.error
## [0, 4) 1 5 4 0 0.90000000
## [4, 6) 0 200 275 0 0.57894737
## [6, 8) 0 100 1988 8 0.05152672
## de 8 a 10 0 0 37 59 0.38541667
RF_importance1 = randomForest::importance(rf_final)[order(randomForest::importance(rf_final)[,1], decreasing = TRUE), ]
plot(rf_final)
legend('topright', colnames(rf$err.rate), col=1:5, fill=1:5)Predição e matriz de confusão - train data
library(caret)
p1 <- predict(rf_final,train1)
head(p1)## 1 2 3 4 5 6
## [6, 8) [6, 8) [6, 8) de 8 a 10 [6, 8) [6, 8)
## Levels: [0, 4) [4, 6) [6, 8) de 8 a 10
head(train1$WR_Grp)## [1] [6, 8) [6, 8) [6, 8) de 8 a 10 [6, 8) [6, 8)
## Levels: [0, 4) [4, 6) [6, 8) de 8 a 10
confusionMatrix(p1, train1$WR_Grp)## Confusion Matrix and Statistics
##
## Reference
## Prediction [0, 4) [4, 6) [6, 8) de 8 a 10
## [0, 4) 10 0 0 0
## [4, 6) 0 474 0 0
## [6, 8) 0 1 2096 0
## de 8 a 10 0 0 0 96
##
## Overall Statistics
##
## Accuracy : 0.9996
## 95% CI : (0.9979, 1)
## No Information Rate : 0.783
## P-Value [Acc > NIR] : < 2.2e-16
##
## Kappa : 0.9989
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: [0, 4) Class: [4, 6) Class: [6, 8)
## Sensitivity 1.000000 0.9979 1.0000
## Specificity 1.000000 1.0000 0.9983
## Pos Pred Value 1.000000 1.0000 0.9995
## Neg Pred Value 1.000000 0.9995 1.0000
## Prevalence 0.003736 0.1774 0.7830
## Detection Rate 0.003736 0.1771 0.7830
## Detection Prevalence 0.003736 0.1771 0.7833
## Balanced Accuracy 1.000000 0.9989 0.9991
## Class: de 8 a 10
## Sensitivity 1.00000
## Specificity 1.00000
## Pos Pred Value 1.00000
## Neg Pred Value 1.00000
## Prevalence 0.03586
## Detection Rate 0.03586
## Detection Prevalence 0.03586
## Balanced Accuracy 1.00000
Predição e matriz de confusão - test data
p2 <- predict(rf_final,test1)
head(p2)## 1 2 3 4 5 6
## [6, 8) [6, 8) [6, 8) de 8 a 10 [6, 8) de 8 a 10
## Levels: [0, 4) [4, 6) [6, 8) de 8 a 10
head(test1$WR_Grp)## [1] [6, 8) [6, 8) [6, 8) [6, 8) [6, 8) [6, 8)
## Levels: [0, 4) [4, 6) [6, 8) de 8 a 10
confusionMatrix(p2, test1$WR_Grp)## Confusion Matrix and Statistics
##
## Reference
## Prediction [0, 4) [4, 6) [6, 8) de 8 a 10
## [0, 4) 0 0 0 0
## [4, 6) 1 79 32 0
## [6, 8) 2 107 841 15
## de 8 a 10 0 0 5 23
##
## Overall Statistics
##
## Accuracy : 0.8534
## 95% CI : (0.8311, 0.8737)
## No Information Rate : 0.7946
## P-Value [Acc > NIR] : 3.007e-07
##
## Kappa : 0.4912
## Mcnemar's Test P-Value : NA
##
## Statistics by Class:
##
## Class: [0, 4) Class: [4, 6) Class: [6, 8)
## Sensitivity 0.000000 0.42473 0.9579
## Specificity 1.000000 0.96409 0.4537
## Pos Pred Value NaN 0.70536 0.8715
## Neg Pred Value 0.997285 0.89225 0.7357
## Prevalence 0.002715 0.16833 0.7946
## Detection Rate 0.000000 0.07149 0.7611
## Detection Prevalence 0.000000 0.10136 0.8733
## Balanced Accuracy 0.500000 0.69441 0.7058
## Class: de 8 a 10
## Sensitivity 0.60526
## Specificity 0.99531
## Pos Pred Value 0.82143
## Neg Pred Value 0.98607
## Prevalence 0.03439
## Detection Rate 0.02081
## Detection Prevalence 0.02534
## Balanced Accuracy 0.80029
Nº de nós nas árvores
hist(treesize(rf_final),
main = "Nº de nós por ávore",
col = "green")Extração de uma única árvore \(Árvore \space n_{tree}=1\)
datatable(getTree(rf_final, 1, labelVar = TRUE),
options = list(searchin = TRUE, pageLength = 5))\(Árvore \space n_{tree}=100\)
datatable(getTree(rf_final, 59, labelVar = TRUE),
options = list(searchin = TRUE, pageLength = 5))remove(INFO)Gráfico de escala multidimensional da matriz de proximidade
MDSplot(rf_final, train1$WR_Grp, pch=20)Plot do modelo rpart - Recursive Partitioning and Regression Trees - personalizando automaticamente a partir do gráfico para o tipo de resposta do modelo.
rp <- rpart::rpart(formula = WR_Grp~.,data=test1)
#rpart::plotcp(rp)
rpart.plot(rp)rpart.plot.version1(rp)O Random Forest apresenta overfiting para os dados de treinamento, como é possível observar na matriz de confusão obtida, e a partir de 60 arvores de regressão no algoritmo, tem-se um ganho minimo ao adicionar novas árvores.
É possível comparar nossos resultados com os resultados da usuária Yueming (disponíveis no site Kaggle clicando AQUI e AQUI), sendo que a partir de seus resultados é possível obter algumas sugestões de melhora do nosso algoritmo, principalmente quanto ao tratamento dos dados (NA e zeros). Em nosso algoritmo, inicialmente, os filmes repetidos não haviam sido considerados, a autora observou 43 observações com repetição e variáveis com excesso de respostas nulas. Quanto as observações repetidas, replicamos o passo sugerido por ela e foram excluídos todos os filmes com repetição.
Em relação as respostas nulas em excesso, ela tratou essas, como NA, entretanto durante as etapas de preparação de dados, descobrimos que a única variável com tal problema era a “movie_facebook_likes” que acabou sendo excluída de nossa análise.
Ainda foram observadas poucos filmes antes de 1980. Tal situação foi enfrentada pela autora com a exclusãos dessas observações raras. Replicamos esse passo também. De uma forma geral, a autora comparou três modelos (KNN, Árvore de Decisão e Random Forest) apenas levando em consideração a acurácia obtida nos modelos.
Aqui, tentamos levar em consideração os valores, principalmente, da matriz de confusão e não somente a acurácia. E podemos concluir que:
table(train2$WR_Grp)##
## [0, 4) [4, 6) [6, 8) de 8 a 10
## 10 475 2096 96
table(test2$WR_Grp)##
## [0, 4) [4, 6) [6, 8) de 8 a 10
## 3 186 878 38
Ainda em análise das classificações corretas vs erradas através da matriz de confusão obtida, com a predição dos dados da base de teste o modelo teve uma acurácia calculada de 85%. E, apesar dos resultados da predição pouco eficiente para todas as categorias, conseguimos obter quais são os fatores/variáveis mais importantes para análise de classificação, dentre as variáveis existentes. Abaixo citamos as top 10 variáveis mais importantes.
RF1 = RF1[1:10,]
hc6_1 <- highchart() %>%
hc_add_series(data = RF1$importance,
type = "bar",
name = "Importância",
showInLegend = FALSE,
tooltip = list(valueDecimals = 2, valuePrefix = "", valueSuffix = ""),
color="orange") %>%
hc_yAxis(title = list(text = "Importância"),
allowDecimals = TRUE, max = 200,
labels = list(format = "{value}")) %>%
hc_xAxis(title = list(text = "Fatores"),
categories = RF1$variables,
tickmarkPlacement = "on",
opposite = FALSE) %>%
hc_title(text = "Importância por fator - Random Forest",
style = list(fontWeight = "bold")) %>%
hc_subtitle(text = paste("")) %>%
hc_tooltip(valueDecimals = 2,
pointFormat = "Importância: {point.y}")%>%
#pointFormat = "Variável: {point.x} <br> Importância: {point.y}")
hc_credits(enabled = TRUE,
text = "Fonte: IMDB/KAGGLE. Elaboração: Ewerson Pimenta.",
style = list(fontSize = "10px")) %>%
hc_exporting(enabled = TRUE, filename = "F6_1-importance-Pimenta")
#hc <- hc %>%
# hc_add_theme(hc_theme_darkunica())
hc6_1Tempo de execução
proc.time() -ptm## user system elapsed
## 47.38 1.59 51.27